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Abstract 

We analyze the Saff man- Taylor viscous fingering problem in rect- 
angular geometry. We investigate the onset of nonlinear effects and 
the basic symmetries of the mode coupling equations, highlighting the 
link between interface asymmetry and viscosity contrast. Symmetry 
breaking occurs through enhanced growth of sub-harmonic perturba- 
tions. Our results explain the absence of finger tip-splitting in the 
early flow stages, and saturation of growth rates compared with the 
predictions of linear stability. 
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1 Introduction 

Pattern formation occurs in many hydrodynamics settings. The Saffman- 
Taylor problem [0, in which two immiscible viscous fluids move in a narrow 
space between two parallel plates (the so-called Hele-Shaw cell), is a widely 
studied example where a fluid-fluid interface evolves @. The initially flat 
interface separating the two fluids can be destabilized by either a pressure 
gradient advancing the less viscous fluid against the more viscous one, or 
by gravity, as a result of a density difference between the fluids. The inter- 
face deforms, and different modes grow and compete dynamically leading to 
undulated patterns. The physics of the first stage of interface instability is 
captured by linear stability analysis ||. After this initial linear behavior, 
the system evolves through a "weakly nonlinear" stage to a complicated late 
stage, characterized by formation of fingers and bubbles, in which nonlinear 
effects dominate. 

The majority of analytical investigations of the dynamics of fingering in- 
stability focus on linear stability analysis or else on selection of steady-state 
patterns [0] . More recent analytical development || addresses the fully non- 
linear time-evolving flow in the small surface tension limit. In contrast, our 

1 



present work develops an analytical approach which links the initial linear 
behavior with the strongly nonlinear advanced stages of the flow. We con- 
sider early stages of the flow, focusing especially on the onset of nonlinear 
effects. We employ an analytical approach known as a mode coupling theory, 
following the method of Haan [|J who studied the Rayleigh- Taylor instabil- 
ity in inertial confinement fusion. We previously applied this method to the 
Saffman- Taylor instability in radial flow geometry Hele-Shaw cells ||. A 
closely related study in the context of solidification was carried out in ref- 
erence ||. Although our approach is quantitatively accurate only at early 
stages of pattern formation, we gain insight into the mechanisms of pattern 
selection and evolution. 

For the rectangular flow geometry, numerical simulations J7|, |5|, |TU 



and experiments |IT, 12[ [13| show that in the nonlinear regime, the viscosity 



contrast A (defined as the difference between the two fluid viscosities divided 
by their sum) plays a major role in the interface evolution. Most noteworthy 
is the role of viscosity contrast in breaking the symmetry of the interface. 
Linear stability analysis does not reveal any interface symmetry breaking in 
the rectangular geometry. In contrast, for the radial geometry the distinction 



between inside and outside of the interface always breaks the symmetry. 

For the rectangular geometry, in the case of zero viscosity contrast (A = 
0), the flow is symmetric about the original flat interface position. On av- 
erage, the fingers of one fluid penetrate the second fluid as much as fingers 
from the second fluid penetrate the first one. For nonzero viscosity contrast 
(A ^ 0), the less viscous fluid penetrates further into the more viscous one. 
The less viscous fluid fingers exhibit stronger length variation than the more 
viscous fluid fingers. For intermediate flow stages, the result is short fingers 
of the more viscous fluid inside the less viscous fluid, and long fingers of the 
less viscous fluid inside the more. In this case, the up-down interface symme- 
try is obviously broken. At late stages the long fingers may pinch off forming 
bubbles. 

In references [0, [8], ||], [K| the mathematical description of the fingering dy- 



namics is written in terms of nonlinear integro-differential equations. These 
equations cannot be solved exactly, and consequently the role of A in asym- 
metry development is primarily revealed by numerically solving these equa- 
tions. A simple analytical description of it remains to be explicitly addressed. 
In addition to the up-down symmetry breaking issue, other possible phe- 



nomena like finger tip-splitting and growth saturation are worth investigat- 
ing. We investigate the saturation of fast-growing modes. By saturation 
we mean a reduction of the exponential rate of growth. Numerical studies 
of Tryggvason and Aref |7], § exhibit a saturation for Saffman- Taylor flow. 
Here, we verify the saturation of growth at third order mode coupling, like 
that seen in the Rayleigh- Taylor instability [f|, [14], [15| . 



Numerical simulations for rectangular geometry flow fT6| , [L7j indicate that 



fingers undergo a type of tip-splitting instability in the late stages of inter- 
face evolution when the dimensionless surface tension parameter B |18| is 
sufficiently small. Experiments |19|, |5D, |2l, 22], performed at very low B, 



observed that the fingers of the less viscous fluid split at their tips while pen- 
etrating the more viscous one. They also observed highly ramified fractal-like 
structures, which hardly resemble a Saffman- Taylor finger. Despite extensive 
numerical and analytical calculations (see, for example, reference 0), consid- 
erable uncertainty surrounds the theoretical understanding of such splitting 
events. In contrast, for the radial flow geometry, finger tip-splitting is com- 
mon experimentally, and it is predicted by a second-order mode coupling 
theory ||. In this paper, we demonstrate that symmetries of the rectan- 



gular flow remove a force driving tip-splitting at second order, explaining 
why it does not commonly occur in early stages of the flow. We investigate 
the possibility of finger tip-splitting at third order, but find that it does not 
occur. 

Section 2 carries out our analytical weakly nonlinear investigation and 
derives a differential equation describing the early nonlinear evolution of the 
interface modes. In section 3, we interpret the results obtained in section 2. 
We identify and analyze the basic symmetries of the mode coupling differ- 
ential equation. The differential equation exhibits the interface asymmetry 
discussed above, and its relation to the value of A. It also indicates the ab- 
sence of finger tip-splitting and the existence of growth saturation at third 
order. We concentrate our attention on the coupling of a small number of 
modes. Section 4 presents our final remarks and perspectives. 

2 The mode coupling differential equation 

Consider two semi-infinite immiscible viscous fluids, flowing in a narrow gap 
of thickness b, in between two parallel plates (see figure 1). We assume that 



5 



b is smaller than any other length scale in the problem, and therefore the 
system is considered to be effectively two-dimensional. Denote the densities 
and viscosities of the lower and upper fluids, respectively as pi, T)\ and p 2 , r] 2 . 
The fluids are assumed to be incompressible, and the flows are assumed to 
be irrotational, except at the interface. Between the two fluids there exists a 
surface tension a. Inject fluid 1 at constant external flow velocity Voo = v^y 
at y = — oo and withdraw fluid 2 at the same velocity at y — +00. We 
describe the system in a frame moving with velocity Voo, so that the interface 
may deform, but it does not displace from y = on the average. In order to 
include the acceleration of gravity g, we tilt the cell so that the y axis lies 
at angle j3 from the vertical direction. During the flow, the interface has a 
perturbed shape described as y = ((x,t) over the range < x < L in the 
comoving frame. 

The early nonlinear evolution of the interface obeys a mode coupling 
equation. We extend Haan's work on the Rayleigh- Taylor problem [|J to the 
case of viscous flow in a Hele-Shaw cell. The analytic model we seek predicts 
the evolution of the interface perturbation amplitude ((x, t). 

We begin by representing the net perturbation ((x,t) in the form of a 
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Fourier expansion 

C(x,t) = ^2Ck(t)exp(ikx), (1) 

k 

where 

(kit) = — J C(x, t) exp(-ikx)dx (2) 

denotes the complex Fourier mode amplitudes. Expansion (^) includes all 
possible modes k, with the exception of k = since we are in a comoving 
frame. The wave vectors are constrained to lie on the x axis, but can be 
either positive or negative. We apply periodic boundary conditions in x 
limiting the wave number k to discrete allowed values 2im/ L, for integer n. 
Experimental realization of the Saffman- Taylor flow in a cylindrical Hele- 
Shaw cell [O] showed that periodic boundary conditions are similar to the 
presence of physical sidewalls. 

In the following paragraphs we use Fourier analysis to gain insight into 
the analytic structure of the flow dynamics. The Fourier approach is of par- 
ticular interest, since despite its simplicity, it allows one to adopt the familiar 
physical ideas of modes and growth rates developed in the linear theory. Here 
we go beyond the level of linear stability analysis and analytically investigate 
the weakly nonlinear stage of the evolution. 



Since we are interested in the early nonlinear behavior of the system, our 
first task is to derive a differential equation for £fc, correct to third order. 
The relevant hydrodynamic equation is Darcy's law [[[]. § 



where Vi = Vi(x,y) and pi = pi(x,y) are, respectively, the velocity and pres- 
sure in fluids i = 1 and 2. Equation @ derives from the Navier-Stokes 
equation by neglecting inertial terms, assuming a parabolic flow profile with 
vanishing velocity at the plates, and by averaging the velocity over the per- 
pendicular direction to the x — y plane. 

— » 

Supplement equation (|3]) with the irrotational flow condition V x vi = 
0. Under this circumstance the velocity is given by = — V0j, where fa 
represents the velocity potential in each of the fluids. Rewrite equation ([|) in 
terms of velocity potentials and then integrate. After dropping an arbitrary 
constant of integration we write Darcy's law for velocity potential 




(3) 



Vi(Pi = {Pi + Pi9V COS 13} + TjiVooy. 



(4) 



Rewrite equation ([|) for each of the fluids and subtract the resulting 
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equations from each other to obtain the jump condition 

\ fx I x I \ of &2 (^2 "Pi) 



12(77! + 



+ C/y (5) 



where 



??2 ~ Vl 

V2 + Vi 



(6) 



is the viscosity contrast and 



, r b 2 (p 2 - pi) g cos p 

U = — ; h AV C 



(7) 



12(77! + V2) 

is a characteristic velocity of the problem. We concentrate attention on 
U > 0, because there is no instability when U < 0. 

The pressure difference at the interface between the two fluids is given 
by § 

fa ~ Pi)\rt = " (Jl + ■ (8) 
The two principal radii of curvature, and describe the interface locally 



The curvature in the x — y plane is p3] 



1 



dx 2 



1 + 



dx 



(9) 



Meanwhile R± ~ 6/2, the radius of curvature in the direction perpendicular 
to the parallel plates, is controlled by the contact angle of the two-fluid in- 



terface at the plates. Typically, one finds that R± is nearly constant [24, 25 



Even though 1/R± 3> l/R\\, the perpendicular curvature does not signifi- 
cantly affect the motion in our problem, since its gradient is nearly zero. 

Combine Darcy's law equation (||), with equations (§) and (J|) for the 
pressure difference and rescale lengths by L and time by L/U, to obtain the 
dimensionless equation of motion 

A (0i|j, =C + 2 | y=C ) - (0!| y=C - (h\y=c) = 2 [C + . (10) 

where 

b 2 a 

~ 12U( Vl + rj 2 )L 2 (U) 
is a dimensionless surface tension coefficient. From now on we work, unless 
otherwise stated, with the dimensionless equation of motion. 

For our weakly nonlinear analysis we are interested in third order con- 
tributions in the perturbation amplitudes. Therefore, all the quantities in 
equation (p~0| ) are evaluated at the perturbed interface y = ((x,t), and not 
at the initial interface position y = as is usually done in linearized surface 
deformation problems. The nonlinear terms arise because of this important 
distinction. 

The incompressibility of the fluids (V • % = 0) means the velocity poten- 
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tials satisfy Laplace's equation V 2 0« = 0. For points far from the interface, 
we expect uniform, steady flow uninfluenced by the interface. Therefore, we 
require that the system evolves with uniform velocity = v^y in the limits 
y — > ±00. Recall that Vi are measured in the comoving frame, and thus 
vanish at y — > ±00. So, 0i and 02 go to constants (that we set to zero) as y 
goes to — 00 and +00, respectively. 

Now define Fourier expansions for the velocity potentials (pi, which must 
obey Laplace's equation V 2 0j = 0, the boundary conditions at y — > ±00, and 
periodic boundary conditions on x. The general velocity potentials obeying 
these requirements are 

0i = J2 ex P(l%) exp(zfcx), (12) 

and 

<p2 = J2 02k (*) exp(-|%) exp(ikx). (13) 

Substitute expansions ( |T2"D and (0) into the equation of motion (]T0|). Keep 
third order terms in the perturbation amplitudes, and then Fourier transform. 
For example, the Fourier transform of the lower fluid velocity potential at 
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the perturbed interface y = ( with wave vector k takes the form 

Mk) = 4>ik(t) + E W\(f>ik'{t)Ck-k' + ~ E (k') 2 <p lk ,(t)( q ( k _ k ,_ q , (14) 

k'^0 1 k',q^0 

where ( k is the Fourier expansion of (, given by equation (||). A similar 
expression for 2 |y=c can be easily obtained. Likewise, the Fourier transform 
of the in-plane curvature with wave vector k, valid up to third order in (, is 

k {l (k) = -k% -f E ( k Tq[k — k' — q]Ck< q Ck-k'- q . (15) 

1 k',q^0 

To close equation (|H]) we need additional relations expressing the velocity 
potentials in terms of the perturbation amplitudes. To find these, consider 
the kinematic boundary condition, which states that the normal components 
of each fluid's velocity at the interface equals the normal velocity of the 



interface itself ||26fl , i.e. 

dt \dxdx) y=( \dy) y=c 



(16) 



Expand equation (|T6"|) to third order in ( and then Fourier transform. Solving 
for (pikif) consistently to third order in £ yields 

<Pik(t) = -j^r + E sgn(kk')Ck'Ck-k' 
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J2 7r\ S 9n{k'q)C, k ,C, q ~k<k-q + ^ TU [ k ~ <? ~ \ J Ck'CqCk-k'-q 
k',q^O l K l k',q^O \ K \ \ A / 

(17) 



and a similar expression for (p2k{t)- Here sgn denotes the sign function. For 
instance, sgn(kk') = 1 if (kk') > and sgn(kk') = — 1 if (kk') < 0. The 
overdot denotes total time derivative. 

Substitute this last expression for 4>ik{t) into equation (|T4D, and again 
keep only cubic terms in the perturbation amplitude (the same procedure 



must be repeated for fluid 2). For k ^ 0, Darcy's law fliCf) becomes the 



differential equation 



(k = \(k)( k + A\k\ Y [1 - sgn(kk')] ( k >(k-k' 

+ \k\\q\sgn(k'q)[l - sgn(kq)]Ck'C q -k'Ck-q 

k',q^0 

h' \h'\\h\ 

k- q- — — -prrr Cfc'CgC 



+ £ k' 

k',q^0 

3 
2 



2 2A; 



qSk — k'—q 



-q- 



k',q^0 



(18) 



where A is the viscosity contrast as defined in equation (|6|), B is the dimen- 
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sionless surface tension coefficient given in equation (|TI ) and 



X(k) = \k\(l - Bk 2 ) (19) 

is the linear growth rate. 

Equation (|18|) is the mode coupling equation of the Saffman- Taylor prob- 
lem for the rectangular geometry Hele-Shaw flow. It gives us the time evolu- 
tion of the perturbation amplitudes (k accurate to third order, conveniently 
written in terms of the two dimensionless parameters A and B. The first 
term on the right-hand side of equation (|i~8D reproduces the linear stability 
analysis [1-3]. The second term, of great importance for understanding the 
interface asymmetry, represents second-order mode coupling. The remaining 
terms depict the third order contributions to the mode coupling equation. 
They lead to saturation of the growth compared with the linear equation of 
motion. Beyond third order, we anticipate that all even terms are multi- 
plied by A, and B enters only into odd terms. In the following section we 
investigate the mode coupling equation (|TJ) in more detail. 
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3 Discussion 

We begin our discussion by analyzing the basic symmetries present in the 
Saffman- Taylor problem and considering how they constrain the form of the 
mode coupling equation ([Tj|). Then we look in detail at features of the 
interface morphology dictated by the first-, second- and third-order terms. 

Several symmetry operations leave the Hele-Shaw cell invariant. Because 
of the periodic boundary conditions introduced in section || the cell is in- 
variant under infinitesimal horizontal translations : (x,y) —>■ (x + e, y). 
Provided the cell is sufficiently long in the vertical direction, we may as- 
sume vertical translational symmetry T* : (x,y) — > (x, y + e). Additional 
symmetries of the cell are the vertical mirror M v : (x,y) (—x,y) and the 
horizontal mirror : (x,y) — ► (x, —y)- 

However, for the Saffman- Taylor problem, the differing fluid densities and 
viscosities, the gravitational acceleration g and external flow break the 
symmetry M^. We define a new symmetry transformation of the Saffman- 
Taylor problem Mh : (x, y, A, B) (x, —y, —A, B). Interchanging viscosities 
r/i and r] 2 reverses the sign of A (see equation (|6])). To understand the trans- 
formation of B defined in equation ( |TTD we must study the transformation 
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of U defined in equation ([/]). Interchanging the densities p\ and P2, and 
simultaneously reversing the direction of gravitational acceleration g leaves 
the first term of equation (|7|) for U invariant. Reversing the sign of A and 
simultaneously reversing the direction of external flow velocity Voo leaves in- 
variant the second term of equation (^) . Invariance of U implies invariance of 
B. Thus Mh swaps the two fluids and reverses the direction of gravitational 
acceleration and external flow at the same time as it reflects the y coordinate. 

The Fourier modes exp(ikx) are basis functions for representations of 
the symmetry group generated by the above operations. By investigating the 
transformation of modes under group generators we can check the symme- 
try properties of the mode coupling equation (|IED, and explain the presence 
and/or absence of various terms. For example, applying horizontal transla- 
tion TJ to a mode of wave vector k multiplies its coefficient Qk by the complex 
phase factor exp(ike). All terms in the mode coupling equation must trans- 
form identically under this translation. Every product of coefficients IlCfcj 

i 

gets multiplied by exp(iJ2kje). Consequently, the wave vectors kj in each 

j 

term must sum up to a common value, k. Inspection of equation (p~8|) verifies 
this rule. 
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Applying vertical translation Ty has no effect (the mode transforms as 
the identity) because ((x) is defined as the interfacial height relative to the 
unperturbed interface. Thus the vertical translation symmetry does not sig- 
nificantly constrain the form of equation (|I~8l) . 

The vertical mirror M v reverses the sign of k. Invariance under this 
operation requires that obeys the same equation as (k- Inspection of 



equation (p~8f ) confirms that all coefficients are even under the reversal of sign 
of all wave vectors. For example, expressions that are first or third order in 
the wave vectors are written as an even power of the wave vectors times an 
absolute value of a wave vector. The even symmetry of the coefficient \(k) 
is explicitly revealed in equation (|T9|). 

The horizontal mirrors (both and M/,) reverse the sign of (k- We 
will discuss the transformation of individual terms in equation ([18|) under 
these mirror operations in subsequent sections as we discuss evolution of the 
interface up to first-, second-, and third-order in our mode coupling theory. 
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3.1 First order - Linear evolution 

The linear stability analysis of the viscous fingering problem has been studied 
since the late 50's |], || . Due to its importance, we briefly review some fea- 
tures of the linear regime. Start with the first order solution to equation (|18|) 

Cf(*) = G(0)exp[A(A:)*]. (20) 

The linear growth rate X(k) (see equation (|19"D ), which leads to exponential 
growth at small amplitudes, is plotted in figure 2 for B = 1. From figure 2 we 
see that, for small wave vectors, perturbations grow in time, deforming the 
interface. At large wave vectors, surface tension stabilizes short wavelength 
deviations. From the linear growth rate ([19]) we can extract two important 
parameters: the critical wave vector (defined by setting X(k) = 0) 

K = -4= (21) 

Vb 



beyond which all modes are linearly stable; and the fastest growing mode 
(defined by setting dX(k)/dk = 0) 

1 (22) 



V3B 

which dominates the initial dynamics of the interface. The fastest growing 
wavelength A* = 2n/k* sets a characteristic length to the problem, giving 



the experimentally observed period of the fingering pattern at initial stages 
of the flow. 

Consider the transformation of the first-order terms in our mode-coupling 
equation fll8|) under horizontal reflections. The coefficient (k reverses sign on 
both sides of the equation. What about X(k)l Since it depends only on 
constant parameters defining the system under study, independent of the 
shape of the interface, it exhibits invariance under M^. As a result the linear 
analysis neither predicts, nor explains, the interface asymmetry observed in 
simulations |7], §, |9|, |10| and experiments JTTJ, [12], [L3j . Although at first-order 



the interface shapes remain (statistically) symmetric under M^, the Saffman- 
Taylor problem does not respect symmetry because of the distinction 
between the lower and upper fluids, and the directions of gravitational ac- 
celeration and external flow. The symmetry operation Mh is the appropriate 
symmetry for the Saffman- Taylor problem. Since B is invariant under M^, 
so is the linear growth rate X(k), confirming the required symmetry. 

It is convenient to rewrite the net perturbation ([!]) in terms of cosine and 
sine modes 

((x, t) — [dk(t) cos(kx) + bk(t) sm(kx)] , (23) 

fc>0 
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where a k = Cfc + C-ft an d °k = i (Ck ~ C-fc) are real-valued. Tryggvason and 
Aref studied the nonlinear behavior of the interface numerically, consid- 
ering the interaction between two cosine waves. Similarly, to illustrate the 
linear evolution of the interface, we perturb the initially flat interface with 
two cosine modes a ks and a kr We consider a dominant fundamental wave 
of wave vector kf = k* = l/y/3 and initial amplitude 0^(0), and a second 
sub-harmonic wave of wave vector k s = kf/2 and relatively weak amplitude 
a ks {0) =0.2a kf {0). 

Exactly two wavelengths of the fundamental mode fit in the Hele-Shaw 
cell, so the fundamental mode is invariant under horizontal translations T^ 2 . 
The sub-harmonic perturbation a ks breaks this translational symmetry of 
the fundamental perturbation a kf by slightly altering the relative lengths 
of the two upward-pointing fingers. Up-down symmetry of the interface is 
weakly broken, by construction, due to our arbitrary choice of sub-harmonic 
O'ks cos (k s x). Had we taken the sub-harmonic bk s sin (k s x) instead, the length 
alternation would have appeared in the downward-pointing fingers. 

Figure 3a depicts the interface given by linear theory. We observe that the 
interface remains nearly up-down symmetric, showing no particular tendency 
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to interface asymmetry development. Each uncoupled mode by itself respects 
up-down symmetry. The translational and reflection symmetry breaking of 
the interface are artifacts of the particular initial condition we select. The 
interfacial deformation grows sufficiently large that quantitative accuracy of 
any perturbative approach is doubtful. We illustrate such large amplitudes 
deliberately, however, because they enhance the visibility of nonlinear effects 
described in the next two subsections. 



3.2 Second order - Role of the asymmetry parameter 

A 

For sufficiently short times we presume that the perturbation series defining 
the mode coupling equation of motion, and its solutions, converge to their 
true forms as successively higher-order terms are incorporated. Truncat- 
ing equation at second order should then result in quantitatively small 
changes in the calculated interface profile ((x,t). Although quantitatively 
small, these changes should incorporate the principal corrections to the lin- 



ear interface evolution discussed above in section |3.1| . Indeed, this section 



explains two important qualitative features of the interfacial evolution: up- 
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down symmetry breaking, and the general absence of finger splitting. 

We consider first the breaking of up-down symmetry. As noted above, 



in section |3.1| , the linear evolution respects reflection symmetry since £k 
reverses sign on both sides of the equation of motion, while A (A;) is unaffected. 
In contrast, the second order term in equation (|l^) breaks Mh symmetry 
because it does not change sign. The second order term is multiplied by 
the viscosity contrast A. Hence the breaking of up-down symmetry depends 
upon the difference between the viscosities of the two fluids. In equation (|18|) 
the role of A in asymmetry development is clearly identified, being revealed 
without solving any complicated nonlinear integro-differential equations [|7|, 
||, |H [Uj . Our mode coupling approach provides a transparent and simple way 
of identifying the intrinsically nonlinear character of the viscosity contrast 
.4. 

Symmetry breaking in the equation of motion translates into symmetry 
breaking in the solution. We now solve equation (|l^) to second order accu- 
racy. Substitute the linear solution (^) into the second-order terms on the 
right hand side of equation (p~8f) , to obtain 
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+ \k\Aj2 [l-sgn(kk')]X(k% k iO)C k MO)ex P [(X(k-k') + X(k'))t] 
+ 0(C* 8 ), (24) 

where A(fc — fc') and A(fc') are the linear growth rates related to the modes 



k — k! and k', respectively. Equation (p4j) is a standard first order differential 



equation [23] with the solution 



C* = Ct +\k\AY, b(0)(k-v(0)G(k, k', t) + 0(( 3 k ), (25) 

where 

A(A0 \ 



g(k,k',t) = [l-sgn(kk') 



X(k - k>) + X(k') - A(Jfe) / 
x (exp[A(A; - + A(fc')]t - exp [A (&)£]) . (26) 

The apparent singularities in Q(k, k',t) at X(k) = X(k') + X(k — k') are can- 
celled by zeros in the numerator of fl26|) , and each term in the sum (equa- 
tion (|25|) ) is regular at these points. At this apparent "resonance point" 
Q(k,k',t) varies as t exp[X(k)t). 

We use the second order solution fl2~5|) to investigate the coupling of a 
small number of modes. As discussed in reference || , even the interaction of 
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two modes can lead to patterns relevant to more complex statistical-fingering 
calculations [0. Our discussion will be simpler if we replace the complex 
Fourier modes ( k with sine and cosine modes as in equation (|2~3"f ). The second 
order equations of motion are 



a k = X(k)a k + Ak^2 [a k >a k+k > + h>b k+ k' , (27) 

fc'>0 

b k = X(k)b k + Ak^2 [^k'bk+k' - b k <a k+k > . (28) 

fc'>0 

Note how the products of sine and cosine amplitudes are arranged to preserve 
the M v symmetry under which the sine functions are odd and the cosines are 
even. Solutions to these equations are similar in form to (|25|). 

Figure 3b illustrates the second-order solution taking the same two-mode 
initial conditions used in figure 3a. Since k s = kf/2, equations (^7|) and (f2lf) 
couple the growth of the sub-harmonic to the amplitude of the fundamental. 
Finger competition, associated with the sub-harmonic mode, is enhanced to 
a degree proportional to the viscosity contrast A. Figure 3b represents the 
extreme case A = 1. Length variation between the upwards fingers is stronger 
than in figure 3a. Notice also that the central upward finger advances more 
strongly into the upper fluid than the downward fingers into the lower fluid. 
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The upward finger at x = advances less strongly. 

Mode coupling not only influences the magnitude of the sub-harmonic, 
but also selects its phase. Without loss of generality we may take a kf > 
and bk f = 0, as we do in the initial conditions for figure 3b. Now compare 
the growth rates of sine and cosine sub-harmonic modes 



Positive viscosity contrast A > increases the growth rate of the cosine sub- 
harmonic a ks , causing increased variability among the lengths of fingers of 
the less viscous lower fluid 1 penetrating downwards into the more viscous 
upper fluid 2. Note that the sign of a ks is dictated by initial conditions, and 
not influenced by mode-coupling. Reversing the sign of ak s has the effect 
of interchanging which of the two upwards pointing fingers will grow at the 
expense of the other. This is tantamount to a horizontal translation by L/2. 
Finally consider patterns with a ka = and b ka ^ 0. Inspecting equation 
we see that A > inhibits the growth of sine modes b ka . Sine modes would 
vary the lengths of fingers of the more viscous upper fluid 2 penetrating into 




(29) 



K k s)h a ~ Ak s b ks a k 



(30) 
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the less viscous lower fluid 1, but their growth is inhibited. 

Reversing the sign of A exactly reverses the above conclusions. Sub- 
harmonic sine modes will be favored over cosine modes. This happens be- 
cause of the Mh invariance of the Saffman- Taylor problem. Reversing the 
signs of both ((x) and A leaves the form of equation (|18D invariant. Conse- 
quently, a randomly chosen up-down symmetric (on average) initial condition 
always evolves into a symmetry-broken interfacial pattern in which the fin- 
gers of the less viscous fluid exhibit variable finger lengths penetrating into 
the more viscous fluid |], |TT], [T^, |13[ . 



Further inspection of the cosine and sine mode coupling equations 
and (^8|) reveals that, while the presence of large wave number modes in- 
fluences the growth of smaller wave number modes, the reverse is not true. 
For example, the growth of any mode k cannot be influenced (up to second 
order) by modes of only smaller wave numbers. In particular, there is no sec- 
ond order term entering the equation of motion for the fundamental mode 
kf in the presence of the sub-harmonic k s . Likewise the presence of kf = k* 
and k s = kf/2 cannot alter the evolution of the harmonic mode kh = 2k f. 
Since kh > k c , the harmonic mode is linearly stable (X(kh) < 0) and will not 
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spontaneously grow. 

This observation yields insight into the absence of finger tip-splitting in 
the rectangular geometry Saffman- Taylor problem, at second order. Splitting 
the tips of the fingers in figure 3b would require the presence of a sizable 
harmonic mode of wave number kh- Even if such modes enter weakly through 
initial conditions or random noise, they quickly die out. 

This is in striking contrast to the radial geometry Saffman- Taylor prob- 
lem ||, where mode coupling drives the growth of harmonic modes. In the 
radial case we may talk about the fundamental cosine mode number rif, 
where rif is an integer counting the number of oscillations around the grow- 
ing perimeter. The influence of a fundamental mode nj on the growth of its 
harmonic mode = 2n f is given by the equation of motion |5[] 

d nh = \{n h )a nh + Cal f , (31) 

where X(rih) denotes the linear growth rate related to the mode (rih) in the 
radial geometry case, and C is negative. Even if \(rih) < 0, the harmonic 
mode can still grow provided the fundamental mode is present, a n . ^ 0. 
In the radial geometry the presence of the fundamental forces growth of 
the harmonic, while in the rectangular geometry the fundamental does not 
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influence the harmonic at second order. Because C < 0, a nh is driven negative, 
the sign that is required to cause fingers to split. 

Consider the behavior of the coefficient of under the radial analogue 
of Mh (here denoted by M™ dm/ ), which interchanges the two fluids and re- 
verses the direction of external flow while reversing the sign of the interfacial 
perturbation (. The coefficient, C, contains both even and odd terms under 
j^radtai _ However, the odd terms are vanishingly small compared to the even 
terms for large unperturbed radii of curvature. We will treat C as effectively 
even under M™ dial in the discussion. 



To transform consistently under Mf l adial , the radial equation of motion fl3T|) 
must contain only terms that are odd. The second order term, however, is 
even and thus breaks M™ dml symmetry. This is expected, because MJ l adial 
is not a symmetry of the radial flow Saffman- Taylor problem. The interface 
is always in-out asymmetric in radial flow because we can always distinguish 
the region of space that lies inside the interface from the region that lies 
outside. For instance, the outward radial motion in which air pushes oil in 
a radial Hele-Shaw cell (divergent flow) is not equivalent to the inward ra- 
dial motion corresponding to withdrawal of oil surrounded by air (convergent 
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flow) [28|| . This broken symmetry allows finger tip-splitting. 

In contrast, Mh is a symmetry of the rectangular Saffman- Taylor prob- 
lem. Thus, even terms in the equation of motion, like the one that causes 
tip-splitting in the radial geometry, are forbidden in the rectangular geom- 
etry. Indeed, it can be shown that the amplitudes of these terms vanish 
when taking the rectangular geometry (large radius) limit of the radial flow 
equations of motion || . At second order the driving force creating finger tip- 
splitting is eliminated by the special symmetry of the rectangular geometry. 

While M h rules out terms like Ca\ f with C even under M h , it does not 
rule out similar terms with odd coefficients such as the viscosity contrast A. 
From our mode coupling equations ( p9|) and fl30|) we know they are absent, 
because small wave number modes do not influence the growth of large wave 
number modes. We lack a symmetry explanation for why this is so. 

Intermediate cases, between radial and rectangular flow geometry, are 
provided by the wedge geometry p8| |29, |30|. In the wedge geometry the 



fluids flow in a Hele-Shaw cell in which the side walls form a wedge with an 
opening angle of 9 (—2n < 9 < 2n), where 9q > (#0 < 0) corresponds 
to a divergent (convergent) flow. Experiments in the wedge geometry [f28 
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with A > observed an increasing sensitivity to finger tip-splitting for larger 
angle 6$ > 0. A related experiment |28|] with A < and convergent flow 



shows that fingers grow, but tip-splitting of inward fingers is inhibited. At 
second order, the specially symmetric case 8q = is unique in the absence of 
a driving force leading to tip splitting. 

3.3 Third order - Onset of saturation and absence of 
tip-splitting 



Now we examine the whole third order mode coupling equation ([Tq) , taking 
into account the contributions coming from the third order terms. Since 
the final expressions are somewhat complicated, we start our discussion by 
considering the evolution of a single mode. The third order mode coupling 
equation reduces to 

"5 



-Bk 2 
2 



Q, (32) 



where we have replaced (k with X(k)(k in terms already of third order on 
the right hand side of equation (|T8D. Since we are interested in the fastest 
growing mode k = k* = 1/ we see from ([32]) that the third order terms 
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lead to a saturation of the growth because the coefficient of £| is negative. 
The exponential growth of the linear instability does not proceed unchecked. 

Figure 3c illustrates the full solution up to third order of equation fll8f) , 
taking the same initial conditions as were used in figures 3a and 3b, and 
including modes a kf / 2 and a*,.. The main effect that is apparent, in compar- 
ison with figure 3b, is significantly diminished amplitude of the fundamental 
mode caused by the saturation effect described by equation fl32|). There is 
a slightly increased amplitude of the sub-harmonic ak s . We also note a slight 
broadening of the dominant central finger, and narrowing of the smaller finger 
at x — 0. 

It is interesting to investigate the possibility of finger tip-splitting at third 
order. Finger tip-splitting is associated to the magnitude of the harmonic 
mode 2k f. It turns out that, at third order the cosines modes a kf and cbk f /2 
force growth of modes a 2 k f and a 3kf / 2 . We consider initial conditions similar 
to those used in figures 3a-c, assuming that modes kf, k s = kf/2 are initially 
present. We study how these two initial modes force growth of modes 2k f and 
3k f/ 2. The influence of the fundamental and sub-harmonic on the growth of 
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the first harmonic kh = 2kf may be expressed by 

3 

< = x ef/a kh - -Bk h k 2 s k f [k f + 2k s ] a kf a 2 ks + 0(af J, (33) 

where 



Xeff = X{k h ) + — 



kfkh 



D (k} + -k 2 h )- 1 



2 



B [ kl + -kl ) - 1 



(34) 

In equation (|33|), the linear growth rate of mode fc^ is increased by the 
presence of kf and k s . Although A(fc/ l ) is negative, opposing the growth 
of the harmonic, the fundamental and sub-harmonic make X e ff less negative. 
Quantitatively, the role of the fundamental dominates the sub-harmonic in 
increasing X e ff- 

It appears at a glance in equation ([£|) that sufficiently large a ks drives 
X e ff positive, permitting growth of the harmonic a kh - It turns out that is 
an artifact of truncating the mode coupling theory. If higher order terms 
in a kf are kept, X e ff will remain negative. Terms involving a kf and a kg 
make A e // less negative but cannot make it go positive. The physical reason 
that prevents positive values for A e // can be understood by considering the 
contour length of the interface. Introducing the harmonic always increases 
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the contour length, although the larger the amplitude of the fundamental, the 
smaller the increase upon introducing the harmonic. Multiplying the contour 
length by the surface tension yields surface energy that favors minimum 
contour length. 

There will be a small amplitude of the harmonic present due to the driving 
term proportional to a,k f a\ a in equation (|33"D , but it is doubtful the amplitude 
will be sufficiently large to split the fingers, because it varies like the second 
power of the sub-harmonic. Therefore, no finger tip-splitting is observed even 
at third order. 

A possible exception to this rule arises from a nonlinear instability caused 
by the coupling of the harmonic mode to itself. Inspecting equation (|32| ) we 



note that high wavenumber modes with k > y 2/5-B (the harmonic mode sat- 
isfies this condition) are unstable for sufficiently large amplitude, because the 
coefficient of is positive. This large-amplitude instability was previously 
noted in the related solidification problem by Dee and Mathur ||. Finger 
tip-splitting in the Saffman- Taylor problem, induced by large amplitude noise 
assisted by a nonlinear instability, has been suggested on the basis of numer- 
ical simulations |16|, [T7| and seen in experiments |19|, [20, |2_T, |22j, performed 
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at very low B. 

4 Concluding remarks 



In this paper we developed a mode coupling theory to investigate the onset 
of nonlinear effects in the viscous fingering problem in a rectangular Hele- 
Shaw cell. From a weakly nonlinear analysis of the system, we derived a 
mode coupling differential equation which describes the evolution of the in- 
terface perturbation amplitudes. The basic symmetries of the mode coupling 
equation are identified and discussed. 

We investigated the relation between the viscosity contrast A and the 
interfacial asymmetry in the Saffman- Taylor problem. Viscosity contrast 
A multiplies symmetry-breaking terms in the mode coupling equation. Our 
analysis explicitly indicates that symmetry breaking occurs through enhanced 
growth of sub-harmonic perturbations. We show that second-order terms 
that drive tip-splitting in the radial flow geometry are prohibited by symme- 
try explaining the general absence of tip-splitting in the rectangular geom- 
etry. A remaining mystery is to explain why small wave number modes do 
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not alter the growth of high wave number modes at second order. 

Our mode coupling analysis shows that modest growth of the harmonic 
may occur at third order by a reduction of the effective linear stability of 
the harmonic mode, coupled with a driving process assisted by finger com- 
petition. We do not expect this process to result in finger tip-splitting under 
normal circumstances. However, self-coupling of the harmonic creates a non- 
linear instability that could split finger tips if noise or other factors drive 
a sufficiently large harmonic amplitude. Finally, we identified the onset of 
saturation effects, which moderates the exponential growth of the linear in- 
stability. 

In a separate work [[31] we extend the present theory to a system in 
which one of the two fluids is a ferrofluid PSI , and a magnetic field is applied 



normal to the Hele-Shaw cell. Interfacial symmetry breaking at late stages 



is very dramatic in this system [32]. We point out, here, that the onset of 



interface symmetry breaking depends on viscosity contrast A, not on the 



applied magnetic field. In reference [31| we show that finger tip-splitting 



arises as a result of the application of an external magnetic field. 
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Figure Captions 



Figure 1: Schematic configuration of the flow in a rectangular Hele-Shaw 
cell. The densities and viscosities of lower and upper fluids are px, r/i and p 2 , 
r/2, respectively. The dashed line represents the unperturbed interface y = 
and the solid undulated curve depicts the perturbed interface y = ((x,t), 
over the range < x < L. The surface tension between the two fluids 
is given by a and the gravitational acceleration, directed in the negative 
vertical direction, is denoted by g. The Hele-Shaw cell of thickness b is tilted 
by an angle (3 from the vertical direction. = v^y represents the uniform 
overall velocity. 



Figure 2: Linear growth X(k) [see equation (|T^) in the text] as a function 
of the wave number k, for B = 1. The critical wave number k c = ±1, and 
the most unstable wave number k* = ±l/-\/3. 

Figure 3: Time evolution of the interface between the fluids, for the case 
of two interacting cosine waves. The dominant mode wave vector kf = k* = 
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1/V3 has initial perturbation amplitude dk f (0) = l. The sub-harmonic mode 
wave vector k s = kf/2 and initial amplitude a,k s (0) = 0.20^(0). For all 
panels the vertical axis scale is the same, the parameter B = 1, and t—0, 
1.5, 3.0 and 4.5; (a) First order (linear) solution; (b) Second order solution 
for A = 1; (c) Third order solution for A = 1. 
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